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Abstract 

The pair-correlation function g{r, t) and its Fourier transform, the structure 
factor S{q,t), are computed during the gelation process of identical spheri- 
cal particles using the diffusion-limited cluster-cluster aggregation model in a 
box. This numerical analysis shows that the time evolution of the character- 
istic cluster size ^ exhibits a crossover close to the gel time tg which depends 
on the volumic fraction c. In this model tg tends to infinity when the box 
size L tends to infinity. For systems of finite size, it is shown numerically 
that, when t < tg, the wave vector qm, at which S{q,t) has a maximum, de- 
creases as S{qm,t)-^/^, where L» is an apparent fractal dimension of clusters, 
as measured from the slope of S{q, t) . The time evolution of the mean num- 
ber of particles per cluster n is also investigated. Our numerical results are 
in qualitative agreement with small angle scattering experiments in several 
systems. 
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I. INTRODUCTION 



In the last years, the kinetics of aggregation has been widely studied |jl|]. It has been 
studied experimentally in a large number of systems like aqueous metal colloids [Q, silicon 
tetramethoxyde and tetraethoxyde [Q, alumino-silcates 0, colloidal silica 0, polystyrene 
|0-|TD[, and oil in water emulsions [TT,T^. Using small-angle scattering techniques it has 



been suggested that some of these systems P,[7|-P|, PT| , P^ grow like diffusion-limited cluster- 
cluster aggregation (DLCA) [|T^|-|T^. Numerically, the gelation process in DLCA has been 
investigated by analyzing the cluster size distribution and the mean cluster size, as a function 
of the time t using both the Smoluchowski equation and the DLCA model in 2 dimensions 



l7| . On the other hand, recently [ll8| , p!9| it has been shown that, when the volumic fraction 



c, or concentration, is greater than a characteristic gel concentration Cg, DLCA model leads 
to a homogeneous gelling network of connected fractal clusters of mean size ^. 

In some systems P,|5|J^^l2| where fractal clusters grow separately until obtaining a system 
of connected fractal clusters at large t, the wave- vector-dependent scattering function S{q, t) 
exhibits a maximum in the wave vector g^-value (which is inversely proportional to ^) 
at a given time t. This maximum can be interpreted as a crossover between the fractal 
regime (intermediate g-values where S{q,t) oc q'^) and the homogeneous regime (small 
g- values). Moreover, some very interesting relations on the S{q) curve (see below equations 
7) have been established on various systems []8|,p!0|. Since scaling relation are known to be 
of wide applications ("universality") it is worth studying them here. Another interesting 
experimental results |^ reveal that the gyration radius Rg (which is proportional to ^) 
saturates close to tg, the time where the system can be considered as a gelling network. 

In this paper, we present some new results obtained from our numerical study of the 
gelation process using diffusion limited cluster-cluster aggregation in a box ||T3HT5|1, where 
a sufficiently large initial concentration has been considered in order to obtain a gelling 
network of connected clusters at the end of the aggregation process. We have calculated 
numerically the mean number of particles per cluster n and their mean size as a function 
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of time t as well as the correlation function g{r,t), and the scattering function S{q,t). We 
show that there exists a crossover at tg such that, for t < tg, C, increases when t increases, 
and for t > tg, the characteristic length ^ remains constant and equal to its value obtained 
near tg. The location of this crossover is pushed toward the low-t values when c increases. 
For c ^ Cg the normalized gel time tg/ttot does not change. Additionally, the maximum 
of the S{q,t) curve, S{qm,t), increases as C,^^' where Dgi is an apparent fractal dimension. 
We report on the discrepancies between the actual fractal dimension D and Dgi for large 
c-values. Finally, we discuss the qualitative agreement between our simulations and some 
experiments. 



II. CONSTRAINTS ON THEORY 



We have modelized the gelation process considering a three dimensional off-lattice exten- 



sion of the original diffusion cluster-cluster aggregation model [0,113 as in ref. [|T^ where we 
have shown that such a model is well adapted to describe the aerogel structure. During the 
aggregation process, we can naturally define a "Monte-Carlo" time t^c which is increased 
by an arbitrary constant amount Stmc at each iteration. Due to the use of formula (3) of 
ref. ||19| to choose the aggregates, the actual "physical" time t should be increased by the 



amount 6t related to 6tmc by: 

5t = (1) 

In our calculation the time unit has been arbitrarily fixed by setting 5tmc = 1- 

Using scaling arguments, it can be established that the whole system remains invariant 
if length, mass and time are conveniently rescaled simultaneously [p0| , pll . Changing the 



average mass n by a factor b, the distance should be rescaled by b^^^ and, using a general 
mean-field argument to evaluate the average time it takes for all clusters to pair up, the 

1 „ (d-2) 

time should be rescaled by 6 -d , giving: 

n (xf^ (2a) 
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with 

For a = —1/D and d = 3 this reasoning gives 7 = 1. In figure 1 we show the time evolution 
of the mean number of particles per cluster n for three different concentrations. When 
c < Cg "h is roughly proportional to t as expected from equation (2). However, for c > Cg 
this dependence is only observed at short times (for t <^tg). At large times (for t 3> tg) if 
one still assumes n ~ P , one should consider a very large exponent 7 ~ 4. Note that in refs. 
20| , |2ll| , equation (2) was obtained under the hypothesis that the entire system remains space 



and time scaling invariant. This hypothesis is no longer valid when a gel is formed since the 
system becomes homogeneous for distances larger than ^ |T8|JT9[| and if one tries to define 



its global fractal dimension it should be taken as the dimension of space. Qualitatively, the 
observed increase of 7 can be understood as a consequence of the increase of the fractal 



dimension. An earlier work by Gonzalez |^2[ suggests the same conclusion. In his paper it 
is shown that when D increases from 1.87 to 2.05, 7 increases from 1.3 to 2.5 when using 
a cluster-cluster aggregation model with a sticking probability p varying from 0.5 to 0.005, 



reaching the chemically-limited cluster-cluster aggregation (CLCA) model PB[, also called 



reaction-limited cluster-cluster aggregation (RLCA) plf. An increasing 7- value with D was 



also found in another CLCA model also recently proposed by Gonzalez . 

The distance that a cluster must travel before colliding with another cluster decreases 
when c increases so that both the gel time tg, when it exists (i. e. for c > c^), and the time at 
the end of the aggregation process ttot should decrease as c increases. This is well observed 
in figure 2 where tg and ttot have been plotted as a function of c in a log-log plot. For c < Cg 
we observe that ttot is roughly inversely proportional to c while, for c > c^, both tg and ttot 



decrease roughly as c"^ ''^ [cg is determined as explained in ref. |T^). The latter result can 
be simply understood if one assumes that when c> Cg the scaling relation (2) holds up to 
t ^ tg where the system is made of clusters of mean size ^ containing particles. The 
concentration being related to ^ by c ~ ^-(3--D) (^^^^ EH), the relation (2) with 7 = 1 



nves: 



D 



tg (X C (3) 

The observed exponent of -1.75 corresponds to D ~ 1.9 a value larger but quite close 
to the expected value of the fractal dimension of DLCA clusters. We point out that in 
our model Cg decreases as L^^^^^^ as showed elsewhere ||19|, and from equation (3) one can 



deduce that tg oo when L oo. From figure 2 we also find that the ratio tg/ttot is almost 
independent on c, for c ^ Cg. This implies that the time taken by smaller clusters to stick 
to the gelling network remains a finite fraction of the total time. 

The two-point correlation function g{r, t) can be calculated by choosing a given path 5r, 
and calculating the number 5n of interparticle distances lying between r and r + 5r, taking 
care of the periodic boundary conditions when investigating regions outside of the box. In 
order to normalize g{r,t) to unity when r tends to infinity, g{r,t) was defined as follows: 

n 1 6n 1 6n , , 

The static structure factor, or scattering function, S{q, t) in a system containing identical 
particles with mean volume fraction c is given by pB| : 



S{q,t) = 1 + - / (^(r,t) - g,)—^ATir'dr (5a) 
vr JO qr 



In theory one should take g^ = 1 and = oo. But, as explained in detail in refs. |[18|J19 
we have calculated g^ from the formula: 

^° - /J- ATTT^dr ^^^^ 

which gives a go value very close to one that insures that S{q,t) when g — * 0. On 
the other hand, was chosen equal to -I, to avoid boundary artifacts due to the periodic 
boundary conditions considered in our problem. 
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III. RESULTS 



Using equation (4), we have calculated g{r, t) at all times t during the aggregation process; 
we have observed that all the g{r,t) curves exhibit a delta peak at r = 1, a discontinuity at 
r = 2 and a minimum at r = ^, that correspond to the typical short and long range features 
of DLCA aggregates [|T^. In figure 3 we have plotted several g{r,t) curves for c = 0.005 at 
four different t emphasizing the region near the minimum in order to show that its location 
(which defines the characteristic size ^ of the aggregates) is shifted to higher r-values when t 
increases. The function g{r,t) was averaged up to 20 independent simulations, and the box 
size L was taken equal to 103. 

Figure 4 shows C, as function of t/ttot for four different values of c. From equation (2), one 
expects that ^ should increase as {t/ttot)° for t < tg. However, for concentrations close to Cg 
this power law is only verified for t -^tg. For t smaller than but close to tg one observes a 
crossover regime with a smaller slope which extends down to smaller times when c is smaller. 
This smaller slope might be due to the interpenetration of clusters which becomes important 
as if: approaches tg. The cluster size ^ saturates for t > tg suggesting that, in this regime, the 
remaining free small clusters diffuse towards the gelling network on which they stick without 
changing the mean size of the large connected clusters. However, during this regime, the 
density of clusters increases and therefore their fractal dimension should become larger. This 
can be checked in figure 5 where we have plotted the fractal dimension D as a function of 
c. The upper curve (open squares) correspond to the "true" fractal dimension, i.e. the one 
estimated in real space from the power law behavior g{r, t) cx r^~^, for 3 < r < ^. 

In figure 6 we show typical S{q,t) curves at different times. These curves exhibit the 
three regimes explained elsewhere fl^, damped oscillations at large g- values (g ^ a^^), a 
power law S{q,t) oc at the intermediate g-values, and a vanishing regime for q <^ 
In this figure we clearly see that is shifted to lower g-values as t increases. This is in 
agreement with the results of figure 4 if one considers that qm is proportional to (as 
found for t = ttot in ref. [0). Additionally, it has been shown ||2^ that the scaled scattering 
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function S{q/qm,t) and the position of the peak are related by: 

S(q/qm, t) = qm{t)-^F{q/q.^) (6) 
where F{q/qm) is a time-independent scahng function. In spinodal decomposition 5 = d = 3, 



while in other growth problems S = D In our case, as showed in inset in figure 6, the 

best collapse is obtained with 6 = 1.9 which roughly corresponds to our previous estimate 
of the fractal dimension. Also, this scaling behavior can be deduced by considering that the 
maximum of the scattering function S{qm, t) is roughly proportional to the mean number of 



particles per cluster |T^: 



S{qrn,t)(xe (7a) 

implying: 

qm oc^ (X S{qm,ty^ (7b) 

In figure 7 we show a log-log plot of |) versus S{qm,t)- For convenience we have 

reported the results for = | by estimating ^ from the minimum of g{r,t) since in such 
estimation the errors are smaller than estimating q^ directly from S{q,t). As expected 
from equation (7b) a quasi linear behavior is found for the four concentrations reported 
in the figure. By performing a linear regression we have estimated the exponents of (7b) 
which have been reported in table I together with the related fractal dimension that we 
call D'. In this table we have also reported the values D^i obtained from the slopes of the 
respective S{q, t) curves at the last stage of the aggregation process as well as the actual 
fractal dimension D obtained from g{r,t). It is shown that, while D increases, D' and Dgi 
decrease when c increases, as already seen on figure 5. Similar discrepancies between D and 
Dgi were observed by Amar et al. |^ in a sub-monolayer molecular beam epitaxy model in 
which the aggregation regime was described by a system of connected fractal clusters. These 
discrepancies between D and D^i might be explained by the fact that for D > 2 the single 
scattering theory considered here is not completely justified Note that, for c> Cg one 



has D = 3 for distances larger than ^, while for c < Cg the final system is a single fractal 



aggregate, with D^i = D a.s showed elsewhere [30 



IV. DISCUSSION 

In this section we would like to discuss our results in comparison with several experi- 
ments. Small-angle neutron-scattering (SANS) experiments in the sol-gel process of alumino- 
silicates have been performed by Pouxviel et al. showing that the radius of gyration Rg 
saturates when t approaches tg (see figure 2 of ref. [^). These authors have considered three 
different composition samples prepared under basic conditions and they obtained a fractal 
dimension D ^ 1.8 as in the DLCA case. Therefore, their results are in agreement with 
our numerical simulations. Another interesting study has been reported in ref. §], where 
small-angle light-scattering experiments have been performed on polystyrene spheres diluted 
in a water-heavy-water mixture; it has been shown that I{qm) scale as qm^, where z is equal 
to -0.58 (close to the exponent reported in table I), at early stage of the growth process, 
and at larger t values it has been found that z becomes close to (-1/3) the typical value for 
spinodal decomposition. They conclude that at early stage of spinodal decomposition the 
growth is dominated by DLCA. Similar kind of analysis has been recently done by Hobbie 



et al. [|lOl concerning SANS experiments in a hydrogen-bonded polymer blend. But in their 
experiments they found D = 2A a, value consistent with diffusion limited particle-cluster 
aggregation ||3ll . 



Other interesting results have been reported by Bibette et al. ||Tl|JT^ . They have per- 
formed small angle light scattering experiments in monodisperse droplets of water emulsions 
in oil. They have suggested that the emulsion growth is like diffusion-limited cluster-cluster 
aggregation, obtaining D 1.8. We have extracted the I{qm,t) points from their curves of 



figure 2 in ref. |Tl|, that we have compared with the points obtained from our simulation 
when using the same c value as in the experiments (c = 0.005). In figure 8 we show the 
good agreement with the experiments, obtaining the same increasing behavior of I{qm,t) as 
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a function of t. 



V. CONCLUSION 

In conclusion, in this work we have shown that numerical simulations of DLCA model in a 
box can account for the growth process of several experimental systems. On the other hand, 
our numerical calculation of S{q, t) during the DLCA process is an original investigation 
that extends previous numerical analysis of this growth process , leading to the 



theoretical explanation of additional interesting properties during the aggregation process 
and confirming some predictions made in some experimental studies. In particular the 
power law given in equation (7b) and the shape of our S{q,t) curves are experimentally 
observed. Similar calculations with the CLCA model are in progress in order to interpret 
some experiments which are more likely explained by such kind of growth process. 
One of us (A. H.) would like to acknowledge support from CONICIT (Venezuela). 
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TABLES 



c exponent (-^) D' Dgi D 

0.005 -0.52 1.92 1.79 1.95 

0.010 -0.54 1.85 1.77 2.02 

0.015 -0.55 1.82 1.76 2.11 

0.030 -0.62 1.61 1.62 2.32 



TABLE I. For four concentrations considered in the simulations, we have reported the respec- 
tive exponents {jjr) from equation (7b) and the related fractal dimensions D', which is compared 
with the apparent fractal dimension Dgi measured from the slope of the S{q, t) curve and the "true" 
fractal dimension D as deduced from g{r,t). 
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FIGURES 

FIG. 1. Log-log plot of n versus ct for L=57.7, and for c = 0.003 (dashed line), c = 0.01 (dotted 
line) and c = 0.05 (solid line). These curves result from averages over 20 simulations. 

FIG. 2. Log-log plot of ttot (black symbols), and tg (open symbols) versus c for L=57.7. The 
solid line corresponds to linear regressions and the vertical arrow show the location of Cg. These 
data result from averages over 20 simulations. 

FIG. 3. Plot of g{r, t) versus r for L = 103 and c = 0.005 at different times t of the aggregation 
process, t/ttot=0-0'i7 (dot-dashed line), t/ttot=0.01 (dashed line), t/ttot=0-57 (dotted line) and 
t/ttot=^ (solid line). These curves result from averages over 10 simulations. 

FIG. 4. Log-log plot of ^ versus t/ttot for L = 103, and for c = 0.005 (open circles), c = 0.01 
(black squares), c = 0.015 (open diamonds) and c = 0.03 (black triangles). The vertical arrows 
indicate the respective location of tg. These data result from averages over 10 simulations. 

FIG. 5. Log-log plot of the "true" fractal dimension D (open symbols) versus c, the open 
squares and open triangles correspond to L = 57.7 and L = 28.85, respectively. The black symbols 
denote the Dgi measured from the slope of the S{q, t) curve for L = 57.7. These curves result from 
averages over 20 simulations. 

FIG. 6. For four different times t/ttot=0.027 (dot-dashed hne), t/ttot=0.01 (dashed line), 
t/ttot=0.57 (dotted line) and t/ttot=^ (solid line) we show the Log-log plot of S{q,t) versus q. 
Inset: Log-log plot of S{q/qm,t)q^. The parameters are the same as in figure 3. 

FIG. 7. Log-log plot of qm versus S{qm,t) for L = 103, and for c = 0.005 (black circles), 
c = 0.01 (open squares), c = 0.015 (open diamonds) and c = 0.03 (black triangles). These curves 
result from averages over 10 simulations. 

FIG. 8. Log-log plot I{qm,t) versus ttot- The open symbols denote the simulation (same pa- 
prameter as in figure 3), and black symbols denote experimental data from ref. [11]. 
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